This invention relates to method and apparatus for imaging a body. In particular, the invention provides a magnetic susceptibility image of an animate or inanimate body. Many imaging techniques exploit some natural phenomenon which varies from tissue to tissue, such as acoustic impedance, nuclear magnetic relaxation, or x-ray attenuation to provide a contrast image of the tissue. Alternatively, some imaging techniques add a substance such as a positron or gamma ray emitter to the body to construct an image of the body by determining the distribution of the substance. Each imaging technique possesses characteristics which result in certain advantages relative to other imaging techniques. For example, the short imaging time of x-ray contrast angiography reduces motion artifacts. In addition, the high resolution of x-ray contrast angiography renders this technique superior to many prior known imaging techniques for high resolution imaging of veins and arteries. However, x-ray contrast angiography is invasive, requires injection of a noxious contrast agent, and results in exposure to ionizing radiation. Thus, it is not typically employed except for patients with severe arterial or venous pathology.
Nuclear Magnetic Imaging (NMR) which is commonly called magnetic resonance imaging (MRI) entails magnetizing a transverse tissue slice with a constant primary magnetic field in a direction perpendicular to the slice, and further magnetizing the slice by applying a gradient in the plane of the slice. A radiofrequency pulse excites selected nuclei of the slice. The excited nuclei relax and emit energy, i.e., radio signals, at frequencies corresponding to local magnetic fields determined by the gradient. A Fourier analysis of the emitted signals provides the signal intensity at each frequency, thereby providing spatial information in one dimension. Repeating the excitation of the nuclei and obtaining the Fourier spectrum of the emitted signals, as the gradient rotates in the plane of the slice, provides a two-dimensional image.
MRI is of primary utility in assessing brain anatomy and pathology. But long NMR relaxation times, a parameter based on how rapidly excited nuclei relax, have prevented NMR from being of utility as a high resolution body imager. The most severe limitation of NMR technology is that for spin echo imaging n, the number of free induction decays (xe2x80x9cFIDsxe2x80x9d), a nuclear radio frequency energy emitting process, must equal the number of lines in the image. A single FID occurs over approximately 0.1 seconds. Not considering the spin/lattice relaxation time, the time for the nuclei to reestablish equilibrium following an RF pulse, which may be seconds, requires an irreducible imaging time of n times 0.1 seconds, which for 512xc3x97512 resolution requires approximately one minute per each two dimensional slice. This represents a multiple of 1500 times longer that the time that would freeze organ movements and avoid image deterioration by motion artifact. For example, to avoid deterioration of cardiac images, the imaging time must not exceed 30 msec. A method for speeding NMR imaging flips the magnetization vector of the nuclei by less than 90 degrees onto the x-y plane, and records less FIDs. Such a method, known as the flash method, can obtain a 128xc3x97128 resolution in approximately 40 seconds. Another technique used to decrease imaging time is to use a field gradient and dynamic phase dispersion, corresponding to rotation of the field gradient, during a single FID to produce imaging times typically of 50 msec. Both methods produce a decreased signal-to-noise ratio (xe2x80x9cSNRxe2x80x9d) relative to spin echo methods. The magnitude of the magnetization vector which links the coil is less for the flash case because the vector is flipped only a few degrees into the xy-plane. The echo-planar technique requires shorter recording times with a concomitant increase in bandwidth and noise. Both methods compensate for decreased SNR by increasing the voxel size with a concomitant decrease in image quality. Physical limitations of these techniques render obtaining high resolution, high contrast vascular images impractical.
Thus, it is an object of the invention to provide high resolution multi-dimensional images of tissue.
It is another object of the invention to provide multi-dimensional magnetic susceptibility images of an object.
It is yet another object of the invention to provide high resolution multi-dimensional images of the cardiopulmonary system.
It is yet another object of the invention to provide a magnetic susceptibility image of a body.
These and other objects of the invention are attained by providing an apparatus for obtaining a multi-dimensional susceptibility image of a body. The apparatus includes a radiation source for magnetizing the body with a magnetic component of a first radiation field. The apparatus also includes a first detector for measuring the magnetic component of the first radiation field in the absence of the body in a volume to be occupied by the body. The apparatus further includes a source for applying a second radiation field to the body, to elicit a third radiation field from the body. A second detector senses this third radiation field, and produces a signal that a reconstruction processor employs to create the magnetic susceptibility image of the body.
One practice of the invention provides a method for determining the distribution of radiation within a magnetized body, emanating from the body in response to an excitation radiation. The method includes the steps of measuring the emanated radiation over a three-dimensional volume by an array of detectors, and uniquely correlating each frequency component of the detected radiation with locations within the body producing that frequency component.
The invention is in part based on the realization that matter having a permeability different from that of free space distorts a magnetic flux applied thereto. This property is called magnetic susceptibility. An object, herein called a phantom, can be considered as a collection of small volume elements, herein referred to as voxels. When a magnetic field is applied to the phantom, each voxel generates a secondary magnetic field at the position of the voxel as well as external to the phantom. The strength of the secondary magnetic field varies according to the strength of the applied field, the magnetic susceptibility of the material within the voxel, and the distance of the external location relative to the voxel. For example, U.S. Pat. No. 5,073,858 of Mills, herein incorporated by reference including the references therein, teaches that the net magnetic flux at a point extrinsic to a phantom to which a magnetic field is applied, is a sum of the applied field and the external contributions from each of the voxels. The ""858 patent further teaches sampling the external flux point by point and employing a reconstruction algorithm, to obtain the magnetic susceptibility of each voxel from the sampled external flux.
Unlike the ""858 patent that relies on a static response from a magnetized body to determine the magnetic susceptibility of the body, a preferred practice of the invention elicits a radiative response from a magnetized body by subjecting the body to a resonant radiation field. One embodiment of the present invention generates a three-dimensional magnetic susceptibility image of an object including a patient placed in a magnetic field from a three-dimensional map of a radio frequency (RF) magnetic field external to the patient, induced by subjecting selected nuclei of the body to a resonant RF field. Application of an RF pulse to the body causes the body to emit the RF magnetic flux external to the body. A Fourier transform of this external flux produces its frequency components (xe2x80x9cLarmor frequenciesxe2x80x9d). Each Larmor frequency determines the magnetic susceptibility of the voxels of the body producing that Larmor frequency. Further, the intensity variation of the external RF field over a three-dimensional volume of space is used to determine the coordinate location of each voxel.
The radiation source for magnetizing a body to be imaged can be a direct current (xe2x80x9cDCxe2x80x9d) magnet, including a superconducting magnet. The radiation sources and amplifiers for applying an RF pulse to the magnetized body are well known in the art, and include, but are not limited to, klystrons, backward wave oscillators, Gunn diodes, and Traveling Wave Tube amplifiers. A preferred embodiment of the invention employs superconducting quantum interference devices (SQUID) as detectors for sensing the external RF field. SQUIDs advantageously allow nulling the applied magnetic flux in order to measure small external fields produced by precessing nuclei. Because the contributions of voxels to an external field at a detector typically drops off as the cube of the distance between the voxel and the detector, some embodiments of the invention employ a magnetizing field whose amplitude increases as the cube of the distance from the detectors. Such a spatial variation of the magnetizing field xe2x80x9clevelsxe2x80x9d the magnitude of the external RF radiation from voxels at different locations relative to the detectors.
One practice of the invention obtains a three-dimensional magnetic susceptibility map of a magnetized body from a three-dimensional map of a secondary magnetic flux produced by the magnetized body, and detected over a three-dimensional volume of space external to the body, herein referred to as the xe2x80x9csample space.xe2x80x9d The extrinsic magnetic flux is sampled at least at the Nyquist rate, i.e., at twice the spatial frequency of the highest frequency of the Fourier transform of the magnetic susceptibility map of the phantom, to allow adequate sampling of spatial variations of the external magnetic flux. This practice of the invention preferably employs a Fourier transform algorithm, described in Fourier Transform Reconstruction Algorithm Section, to form the magnetic susceptibility map of the object.
One embodiment of the present invention employs nuclear magnetic resonance (NMR) to induce a magnetized body to emit an external radiation having a magnetic field component. In particular, application of an RF pulse, resonant with selected nuclei of a magnetized body, can polarize the nuclei through rotation of their magnetic moments. The polarized nuclei within a voxel precess about the local magnetic field in the voxel at a Larmor frequency determined in part by the applied magnetic field at position of the voxel and the susceptibility of the voxel. The superposition of external RF fields produced by all the voxels of the body creates the total external RF field. Thus, the external RF field contains frequency components corresponding to precession frequencies of the nuclei in different voxels of the body.
One practice of the invention detects the external RF field in the near field region where the distance of a detector sensing radiation from a voxel at a distancer from the detector is much smaller than the wavelength xcex of the radiation emitted by the voxel, i.e., r less than  less than xcex (or kr less than  less than 1). The near fields are quasi-stationary, that is they oscillate harmonically as exe2x88x92ixcfx89t, but are otherwise static in character. Thus, the transverse RF magnetic field of each voxel is that of a dipole. In one embodiment, an array of miniature RF antennas sample the external RF field over a three-dimensional volume of space that can be either above or below the object to be imaged. The frequency components of the detected RF signals determine the magnetic susceptibility of voxels that give rise to the RF signals, and the location of each voxel is determined through measurements of the spatial variations of the intensity of the external RF field at a given frequency. Thus, the frequency components of the external RF radiation, and the intensity variations of the external RF radiation provide the necessary information for providing a magnetic susceptibility map of the magnetized phantom, such as a human body. Such a susceptibility map can be employed to obtain anatomical images of a human body based on selected physiological parameters.
Another practice of the invention employs paramagnetic and/or diamagnetic substances present in a body to be imaged to cause variations of local magnetic fields in different parts of the body, thereby providing a susceptibility image of the body. In particular, an excitation field, such a magnetic field, polarizes the paramagnetic and/or diamagnetic substances of the body. The polarized substances contribute to the local magnetic field at the positions of the voxels comprising the body. The amount of this contribution varies from one voxel to another depending on the variation of the concentration of the substances throughout the body. An excitation of selected nuclei of the voxels cause the nuclei to provide an external RF field through precession about the local magnetic fields. The external RF field includes a range of frequencies relating to the local magnetic field at position of each voxel. Thus, an analysis of the frequencies lead to information regarding the distribution of the paramagnetic or diamagnetic substances throughout the body.
All matter has a permeability different than that of free space that distorts an applied magnetic flux. This property is called magnetic susceptibility. An object, herein called a phantom, can be considered as a collection of small volume elements (hereafter called xe2x80x9cvoxelsxe2x80x9d). When a magnetic flux is applied to the phantom, each voxel generates a secondary magnetic flux at the position of the voxel as well as external to the phantom. The strength of the secondary magnetic flux varies according to the strength of the applied flux, the magnetic susceptibility of the material within the voxel, and the distance of the external location relative to the voxel. In one embodiment described in Mills [1] which is herein incorporated by reference including the references given therein in the Reference Section, the net magnetic flux at a point extrinsic to the phantom is a sum of the applied flux and the contributions of each of the voxels within the object. This flux is point sampled over a three dimensional space, and the magnetic susceptibility of each voxel is solved by a reconstruction algorithm.
In an embodiment of the present invention, the three-dimensional magnetic susceptibility map of an object including a patient placed in a magnetic field is generated from a three-dimensional map of the transverse resonant radio frequency (RF) magnetic flux external to the patient. The magnetic susceptibility of each voxel is determined from the shift of the Larmor frequency due to the presence of the voxel in the magnetizing field. The intensity variation of the transverse RF field over space is used to determine the coordinate location of each voxel. The RF field is the near field which is a dipole that serves as a basis element to form a unique reconstruction. The geometric system function corresponding to a dipole which determines the spatial intensity variations of the RF field is a band-pass for kxcfx81=kz. In the limit, each volume element is reconstructed independently in parallel with all other volume elements such that the scan time is no greater than the nuclear free induction decay (FID) time.
The magnetic moment, mZ, of each voxel within the phantom is given by the product of its volume magnetic susceptibility, "khgr", the magnetizing flux oriented along the z-axis, BZ, and the volume of the voxel, xcex94V.
mZ="khgr"Bxcex94V,xe2x80x83xe2x80x83(1)
The magnetic moment of each voxel is a magnetic dipole. And the phantom can be considered to be a three-dimensional array of magnetic dipoles. At any point extrinsic to the phantom, the z-component of the secondary flux, Bxe2x80x2, from any single voxel is                               B          xe2x80x2                =                                            m              z                        ⁢            2            ⁢                          z              2                                -                      x            2                    -                                    y              2                                                      (                                                      x                    2                                    +                                      y                    2                                    +                                      z                    2                                                  )                                            5                /                2                                                                        (        2        )            
where x, y, and z are the distances from the center of the voxel to the sampling point. It is shown below that no geometric distribution of magnetic dipoles can give rise to Eq. (2). Therefore, the flux of each magnetic dipole (voxel contribution) forms a basis set for the flux of the array of dipoles which comprise the magnetic susceptibility map of the phantom.
Eq. (2) is a system function which gives the magnetic flux output in response to a magnetic dipole input at the origin. The phantom is an array of spatially advanced and delayed dipoles weighted according to the magnetic susceptibility of each voxel; this is the input function. The secondary flux is the superposition of spatially advanced and delayed flux, according to Eq. (2); this is the output function. Thus, the response of space to a magnetized phantom is given by the convolution of Eq. (2) with the series of weighted, spatially advanced and delayed dipoles representing the susceptibility map of the phantom.
In Fourier space, the output function is the product of the Fourier transform (FT) of the system function and the FT of the input function. Thus, the system function filters the input function. The output function is the flux over all space. However, virtually all of the spectrum (information needed to reconstruct the magnetic susceptibility map) of the phantom exists in the space outside of the phantom because the system function is essentially a band-pass filter. This can be appreciated by considering the FT, H[kxcfx81, kZ], of Eq. (2):                               H          ⁡                      [                                          k                ρ                            ,                              k                z                                      ]                          =                                            4              ⁢              π              ⁢                              xe2x80x83                            ⁢                              k                ρ                2                                                                    k                z                2                            +                              k                ρ                2                                              =                                    4              ⁢              π                                      1              +                                                k                  z                  2                                                  k                  ρ                  2                                                                                        (        3        )            
where kxcfx81 is the spatial frequency in the xy-plane or kxcfx81-plane and kZ is the spatial frequency along the z-axis. H[kxcfx81, kZ] is a constant for kxcfx81 and kZ essentially equal as demonstrated graphically in FIG. 1c. 
A magnetic field with lines in the direction of the z-axis applied to an object comprised of magnetically susceptible material magnetizes the material so that a secondary field superposes the applied field as shown in FIG. 9. The secondary field outside of the object (phantom) and detected at a detector 301 is that of a series of magnetic dipoles centered on volume elements 302 of the magnetized material. In Cartesian coordinates, the secondary magnetic flux, Bxe2x80x2, at the point (x,y,z) due to a magnetic dipole having a magnetic dipole moment mz at the position (x0,y0,z0) is                               B          xe2x80x2                =                              μ            0                    ⁢                                                    m                z                            ⁡                              (                                                      2                    ⁢                                                                  (                                                  z                          -                                                      z                            0                                                                          )                                            2                                                        -                                                            (                                              x                        -                                                  x                          0                                                                    )                                        2                                    -                                                            (                                              y                        -                                                  y                          0                                                                    )                                        2                                                  )                                                                    [                                                                            (                                              x                        -                                                  x                          0                                                                    )                                        2                                    +                                                            (                                              y                        -                                                  y                          0                                                                    )                                        2                                    +                                                            (                                              z                        -                                                  z                          0                                                                    )                                        2                                                  ]                                            5                /                2                                              ⁢                      i            z                                              (        4        )                                          B          xe2x80x2                =                                                            (                                                      2                    ⁢                                          z                      2                                                        -                                      x                    2                                    -                                      y                    2                                                  )                                                              [                                                            x                      2                                        +                                          y                      2                                        +                                          z                      2                                                        ]                                                  5                  /                  2                                                      ⊗                          m              z                                ⁢                      δ            ⁡                          (                                                x                  -                                      x                    0                                                  ,                                  y                  -                                      y                    0                                                  ,                                  z                  -                                      z                    0                                                              )                                ⁢                      i            z                                              (        5        )            
where iZ is the unit vector along the z-axis. The field is the convolution of the system function h(x,y,z) or h(xcfx81,"PHgr",z) (the left-handed part of Eq. (5)), with the delta function (the right-hand part of Eq. (5)), at the position (x0,y0,z0). A very important theorem of Fourier analysis states that the Fourier transform of a convolution is the product of the individual Fourier transforms [2]. The Fourier transform of the system function, h(x,y,z) or h(xcfx81,"PHgr",z) is given in APPENDIX V.
The z-component of a magnetic dipole oriented in the z-direction has the system function, h(x,y,z), which has the Fourier transform, H[kx,ky,kz], which is shown in FIG. 1c.                                                                         H                ⁡                                  [                                                            k                      x                                        ,                                          k                      y                                        ,                                          k                      z                                                        ]                                            =                              xe2x80x83                            ⁢                                                4                  ⁢                                      π                    ⁡                                          [                                                                        k                          x                          2                                                +                                                  k                          y                          2                                                                    ]                                                                                        [                                                            k                      x                      2                                        +                                          k                      y                      2                                        +                                          k                      z                      2                                                        ]                                                                                                        xe2x80x83                            ⁢                              (                6                )                                                                                        =                              xe2x80x83                            ⁢                                                H                  ⁡                                      [                                                                  k                        ρ                                            ,                                              k                        z                                                              ]                                                  =                                                                            4                      ⁢                      π                      ⁢                                              xe2x80x83                                            ⁢                                              k                        ρ                        2                                                                                                            k                        z                        2                                            +                                              k                        ρ                        2                                                                              =                                                            4                      ⁢                      π                                                              1                      +                                                                        k                          z                          2                                                                          k                          ρ                          2                                                                                                                                                                                            xe2x80x83                            ⁢                              (                7                )                                                                        xe2x80x83            
The output function, the secondary magnetic field, is the convolution of the system function, h(x,y,z)xe2x80x94the geometric transfer function for the z-component of a z-oriented magnetic dipole with the input functionxe2x80x94a periodic array of delta functions each at the position of a magnetic dipole corresponding to a magnetized volume element.                                           (                                          2                ⁢                                  z                  2                                            -                              x                2                            -                              y                2                                      )                                              [                                                x                  2                                +                                  y                  2                                +                                  z                  2                                            ]                                      5              /              2                                      ⊗                              ∑                          n              =                              -                ∞                                      ∞                    ⁢                      xe2x80x83                    ⁢                                    m              z                        ⁢                          δ              ⁡                              (                                                      x                    -                                          nx                      0                                                        ,                                      y                    -                                          ny                      0                                                        ,                                      z                    -                                          nz                      0                                                                      )                                                                        (        8        )            
The Fourier transform of a periodic array of delta functions (the right-hand side of Eq. (8)) is also a periodic array of delta functions in k-space:                               1                                    x              0                        ⁢                          y              0                        ⁢                          z              0                                      ⁢                              ∑                          n              =                              -                ∞                                      ∞                    ⁢                      xe2x80x83                    ⁢                                    m              z                        ⁢                          δ              ⁡                              (                                                                            k                      x                                        -                                          n                                              x                        0                                                                              ,                                                            k                      y                                        -                                          n                                              y                        0                                                                              ,                                                            k                      z                                        -                                          n                                              z                        0                                                                                            )                                                                        (        9        )            
By the Fourier Theorem, the Fourier transform of the spatial output function, Eq. (8), is the product of the Fourier transform of the system function given by Eq. (7), and the Fourier transform of the input function given by Eq. (9).                                           4            ⁢            π                                1            +                                          k                z                2                                            k                ρ                2                                                    ⁢                  1                                    x              0                        ⁢                          y              0                        ⁢                          z              0                                      ⁢                              ∑                          n              =                              -                ∞                                      ∞                    ⁢                      xe2x80x83                    ⁢                                    m              z                        ⁢                          δ              ⁡                              (                                                                            k                      x                                        -                                          n                                              x                        0                                                                              ,                                                            k                      y                                        -                                          n                                              y                        0                                                                              ,                                                            k                      z                                        -                                          n                                              z                        0                                                                                            )                                                                        (        10        )            
In the special case that
kxcfx81=kzxe2x80x83xe2x80x83(11)
the Fourier transform of the system function (the left-hand side of Eq. (10)) is given by
H=4xcfx80xe2x80x83xe2x80x83(12)
Thus, the Fourier transform of the system function band-passes the Fourier transform of the input function. Both the input function (the right-hand part of Eq. (8)) and its Fourier transform (the right-hand part of Eq. (10)) are a periodic array of delta functions. No frequencies of the Fourier transform of the input function are attenuated; thus, no information is lost in the case where Eq. (11) holds. Thus, the resolution of the reconstructed magnetic susceptibility map is limited by the spatial sampling rate of the secondary magnetic field according to the Nyquist Sampling Theorem.
Superconducting Quantum Interference Devices (SQUIDs) allow the user to null the ambient or applied magnetic flux and to measure very small secondary magnetic fluxes (10xe2x88x9215 Tesla). The flux contribution from a voxel of human tissue of millimeter dimensions is typically 10xe2x88x929 Tesla. For example, representative parameters of Eq. (1) are "khgr"=10xe2x88x924, BZ=105G, and xcex94V=10xe2x88x923 cm3 which results in mz=10xe2x88x922 Gcm3. Substitution of this voxel magnetic moment into Eq. (2) with x=0; y=0; z=10 cm results in B=10xe2x88x925 G=1 nT. Eq. (2) shows that the field at the detector drops off as the cube of the distance from the voxel of tissue. By applying a voxel magnetizing field that increases as the cube of the distance from the detector to the voxel, the field at the detector becomes independent of the distance from the voxel. In this case, the field is given by
Bz=B0[a2+zn2]3/2xe2x80x83xe2x80x83(13)
where a0 and B0 are constants and zn is the distance from the detector to the voxel. Thus, a magnetizing field that increases as the distance cubed from the detector xe2x80x9clevelsxe2x80x9d the magnitude of signals from voxels at different locations relative to the detector. A three-dimensional magnetic susceptibility map can be generated from a three-dimensional map of the secondary magnetic flux of a magnetized phantomxe2x80x94such as a human (the corresponding space of detection is hereafter called the xe2x80x9csample spacexe2x80x9d). The magnetic susceptibility map is reconstituted using a Fourier transform algorithm. The algorithm is based on a closed-form solution of the inverse problemxe2x80x94solving the spatial distribution of magnetic dipoles from the extrinsic secondary flux from an array of magnetized voxels. The extrinsic magnetic flux is sampled at the Nyquist rate (twice the spatial frequency of the highest frequency of the Fourier transform of the magnetic susceptibility map of the phantom).
In an embodiment of the present invention, nuclear magnetic resonance (NMR) is the means to measure the secondary magnetic field to provide the input to the magnetic susceptibility reconstruction algorithm. In this case, the measured secondary (RF) field is transverse to the magnetic flux including the local contribution due to the magnetic susceptibility of the voxel. The RF field is detected in the near zone where r less than  less than xcex (or kr less than  less than 1), and the near fields are quasi-stationary, oscillating harmonically as exe2x88x92ixcfx89t, but otherwise static in character. Thus, the transverse RF magnetic field of each voxel is that of a dipole, the maximum amplitude of which is given by Eq. (2) wherein the Larmor frequency of each voxel is shifted due to its magnetic susceptibility, and mz, the magnetic moment along the z-axis, of Eq. (2) corresponds to the bulk magnetization M of each voxel. In terms of the coordinates of Eq. (2), an array of miniature RF antennas point samples the maximum dipole component of the RF signal over the sample space such as the half space above (below) the object to be imaged wherein each RF signal is frequency-shifted by the perpendicularly oriented magnetic susceptibility moment of each voxel. A three-dimensional magnetic susceptibility map is generated from a three-dimensional map of the secondary (RF) magnetic flux of a magnetized phantomxe2x80x94such as a human. The magnetic susceptibility of each voxel(s) is determined from the shift of the Larmor frequency due to the presence of the voxel(s) in the magnetizing field. The measurements of the spatial variations of the intensity of the transverse RF field of a given frequency is used to determine the coordinate location of each voxel(s). In one embodiment, the magnetic susceptibility map is reconstructed using a Fourier transform algorithm. The algorithm is based on a closed-form solution of the inverse problemxe2x80x94solving the spatial location(s) of a magnetic dipole(s) of known magnetic moment (via the Larmor frequency) from the spatial variations of the extrinsic transverse secondary (RF) flux from the magnetized voxel(s). The extrinsic (RF) magnetic flux is sampled at the Nyquist rate (twice the spatial frequency of the highest frequency of the Fourier transform of the magnetic susceptibility map of the phantom). In the limit, each volume element is reconstructed independently in parallel with all other volume elements such that the scan time is no greater than the nuclear free induction decay (FID) time.
The NMR type ReMSI scan performed on the object to be imaged including a human comprises the following steps:
The primary or magnetizing field is first determined over the voxels of the imaged volume (hereafter called xe2x80x9cimage spacexe2x80x9d) in the absence of the object to be imaged.
The magnetic moments of nuclei including protons of the object to be imaged that are aligned by the primary field are further aligned by a radio frequency (RF) pulse or series of pulses.
The strength and duration of the rotating H1 (RF) field that is resonant with the protons of the magnetized volume and is oriented perpendicularly to the direction of the magnetizing field is applied such that the final precession angle of the magnetization is 90xc2x0 (xcfx86=90xc2x0) such that the RF dipole is transverse to the primary magnetizing field and perpendicular to the RF magnetic field detector.
The nuclei precess and emit RF energy at the frequency determined by the total magnetic flux at the voxel which is the superposition of the primary (applied flux) and the secondary (magnetic susceptibility flux).
The precessing nuclei undergo a free induction decay as they emit energy.
The magnetic intensity is time harmonic according to the magnetic flux at each voxel. The time dependent signals are Fourier transformed to give the spectrum. The magnetic moment corresponding to each frequency is determined (Eq. (21)).
The transverse RF field is a dipole. The maximum intensity is determined for each frequency at each detector for at least one synchronized time point or short time interval to form a matrix for each frequency of the intensity as a function detector position in sample space. The measurements of the spatial variations of the transverse RF field of a given frequency is used to determine the coordinate location of each voxel(s). Thus, each matrix comprises the intensity variation over the sample space of the RF field of the bulk magnetization M of each voxel having the corresponding magnetic moment. To improve the signal to noise ratio of the image, this procedure may be repeated for multiple time points with reconstruction and averaging of the superimposed images, or the data may be averaged before reconstruction.
The Fourier transform algorithm is performed on each frequency component over the detector array to map each bulk magnetization M corresponding to a magnetic moment to a spatial location or locations over the image space. The later case applies if the same magnetic moment is at more than one position in the image. The three-dimensional magnetic susceptibility map (the input function) is given as follows. With respect to the coordinate system of Eq. (2), (x, y, and z are the Cartesian coordinates, mz, the magnetic moment along the z-axis, of Eq. (2) corresponds to the bulk magnetization M of each voxel, and B is the magnetic flux due to the magnetic moment shown in FIG. 9; the relationship to the NMR coordinate system is given in the Reconstruction Algorithm Section) the origin of the coordinate system, (0,0,0), is the center of the upper edge of the phantom. The phantom occupies the space below the plane x, y, z=0 (zxe2x89xa60) in the phantom space), and the sampling points lie above the plane (z greater than 0 in the sampling space). The magnetic flux in the sampling space is given by multiplying the convolution of the input function with the system function by the unitary function (one for zxe2x89xa70 and zero for z less than 0). The input function can be solved in closed-form via the following operations:
1. Measure the magnetic flux at discrete points in the sampling space. Each point is designated (x, y, z, flux) and each flux value is an element in matrix A.
2. Discrete Fourier transform matrix A to obtain matrix B.
3. Multiply each element (flux value) of matrix B by the corresponding inverse (reciprocal) value of the Fourier transform of the system function, Eq. (2), evaluated at the same frequency as the element of the matrix A. This is matrix C.
4. Generate matrix D by taking the discrete inverse Fourier transform of matrix C.
5. Multiply each element of matrix D by the distance squared along the z-axis to which the element corresponds to generate the bulk magnetization M map. (This corrects the limitation of the sample space to zxe2x89xa70).
6. In one embodiment, the voxels with a finite bulk magnetization M above a certain threshold or at an edge in the map are assigned the magnetic susceptibility (magnetic moment) determined from the shift of the Larmor frequency due to the presence of the voxel(s) in the uniform (nonuniform) magnetizing field. The other voxels are assigned a zero value. This procedure is repeated for all Larmor frequencies. In the limit with sufficient Larmor frequency resolution, each volume element is reconstructed independently in parallel with all other volume elements such that the scan time is no greater than the nuclear free induction decay (FID) time.
7. In the case that the primary field is nonuniform, the magnetic moment map is determined, then the magnetic susceptibility map of each magnetic moment is given by dividing the magnetic moment map by the primary magnetic field map on a voxel by voxel basis and by subtracting one from each term (Eq. (23)). The total magnetic susceptibility map is the superposition of the separate maps for each magnetic moment corresponding to a unique Larmor frequency which is plotted and displayed.